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Abstract: Precision tests of gravity can be used to constrain the properties of hypothetical 
very hght scalar fields, but these tests depend crucially on how macroscopic astrophysical 
objects couple to the new scalar field. We study the equations of stellar structure using 
scalar-tensor gravity, with the goal of seeing how stellar properties depend on assumptions 
made about the scalar coupling at a microscopic level. In order to make the study relatively 
easy for different assumptions about microscopic couplings, we develop quasi- analytic ap- 
proximate methods for solving the stellar-structure equations rather than simply integrating 
them numerically. (The approximation involved assumes the dimensionless scalar coupling 
at the stellar center is weak, and we compare our results with numerical integration in or- 
der to establish its domain of validity.) We illustrate these methods by applying them to 
Brans-Dicke scalars, and their generalization in which the scalar- matter coupling slowly runs 
- or 'walks' - as a function of the scalar field: a{(j)) ~ a<j -|- hs4>- (Such couplings can arise in 
extra-dimensional applications, for instance.) The four observable parameters that charac- 
terize the fields external to a spherically symmetric star are the stellar radius, R, mass, M, 
scalar 'charge', Q, and the scalar's asymptotic value, 0oo- These are subject to two relations 
because of the matching to the interior solution, generalizing the usual mass-radius, M{R), 
relation of General Relativity. S ince is common to different stars in a given region (such 
as a binary pulsar), all quantities can be computed locally in terms of the stellar masses. 
We identify how these relations depend on the microscopic scalar couplings, agreeing with 
earlier workers when comparisons are possible. Explicit analytical solutions are obtained for 
the instructive toy model of constant-density stars, whose properties we compare to more 
realistic equations of state for neutron star models. 
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1. Introduction 

One way in which candidate quantum theories of gravity can differ from one another is the 
spectrum of bosons they predict at very low energies. Such fields, if they exist and are 
sufficiently light, can mediate long-range forces that can observably compete with gravity. 

Scalar-tensor theories of gravity are among those that can arise in this way, with light 
scalar fields in addition to the usual metric tensor [|l|, |2|, ^, ^ ^, A variety of scalars 
commonly arise in fundamental theories, although it is very unusual for them to be light 
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enough to mediate forces over macroscopic distances. They are rarely this hght because 
quantum corrections famously tend to give scalars masses, even if they would have been 
massless at the purely classical level. But in some circumstances symmetries can protect 
against masses, such as if the scalar is a pseudo-Goldstone boson Q for a spontaneously 
broken approximate symmetry, or part of the low-energy limit of an extra-dimensional model 

i- 

A good deal of effort has been invested in comparing the predictions of scalar-tensor 
theories with observations in various astrophysical systems ||9|, 11, 12, 14, |15|, 17], 



in order to improve the constraints on their existence (or to discover their presence). Binary 
pulsars are particularly useful for this purpose, since the precision of their timing allows 
accurate measurements of the relativistic gravitational effects that are generated by the strong 
gravitational fields present. 

Central to these tests is an understanding of how the stellar properties depend on the 
microscopic couplings of the scalar field to matter. Yet these stellar properties can sometimes 
be remarkable in scalar-tensor theories in the presence of relativistic sources. In particular 
a phenomenon called spontaneous scalarization can occur, in which the star above a critical 
baryonic mass can locally support a scalar field even if the scalar-matter coupling vanishes 



asymptotically far from the star |18]. This allows scalar-tensor gravity to deviate strongly 
from General Relativity (GR) near the star while still passing weak-field solar system tests. 

Since different fundamental theories predict scalars with different microscopic couplings, 
it is useful to be able to survey how stellar properties depend on these couplings. For this rea- 
son in this paper we re-examine the equilibrium conditions in a star in scalar-tensor models as 
a function of scalar couplings. In particular we do so working as far as possible within analytic 
approximation schemes, since these more easily allow the results to be varied for different kinds 
of scalar properties. We find semi-analytic progress is possible using a weak-coupling approx- 
imation for the scalar field near the stellar center. By comparing with numerical integrations 
of the equations of stellar structure we find the domain of validity of these approximations, 
which allow us to understand stellar properties fairly well. 

As a preliminary to studying the relativistic limit in more detail, we apply our analysis to 
the special case of an incompressible stellar fluid, for which the energy density is approximately 
constant. (This can be done consistent with conservation laws only if the pressure is not also 
regarded as being a function of the energy density.) To leading order in the weak-coupling 
expansion, the profiles of physical variables of constant-density stars can be described in 
closed form in terms of elementary functions, dilogarithms, and Heun functions. 

Although not physically realistic, solutions with this equation of state can give insight 
into general features of relativistic stars, in particular near the maximum mass that can be 
supported by gravitational forces. (For GR the maximum mass found under the assumption 
of constant density gives an upper bound on the maximum mass that would be found with 



more realistic equations of state [19[.) We compare the predictions of constant density in 
scalar-tensor theory with those of several representative equations of states for neutron stars. 
The rest of the paper is organized as follows. In §2 we define the scalar-tensor theories 
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of interest, and in particular the parameters describing the couplings of the light scalar to 
ordinary matter. Following much of the literature we specialize to the case where the scalar- 
matter couplings are only weakly dependent on the scalar fields themselves — what we call 
quasi-Brans/Dicke (qBD) models — as well as the constraints on their couplings that are 
inferred from solar system tests of gravity. §3 then derives the equations describing hydro- 
static equilibrium for static and spherically-symmetric stars in scalar-tensor gravity, as well 
as the matching formulae that relate the interior and exterior geometries. Some of those 
general properties that can be extracted without solving them explicitly are also discussed, 
such as whether the pressure need be monotonic; the kinds of relations they imply amongst 
quantities visible to exterior observers; and the non-relativistic limit. Next, §4 provides a 
discussion of the perturbative solutions to these equations in the weak-coupling limit, in both 
the relativistic and non-relativistic cases. The important distinction between perturbations 
in the strength of the coupling measured at the stellar center, vs its strength at infinity, first 
arises in this section. Finally, in section §5 the perturbed equations are solved explicitly to 
leading nontrivial order for the special case of incompressible stars, with the results compared 
to more realistic equations of state for neutron star models. 

2. Single-field scalar-tensor models 

We start by defining the field equations of the scalar-tensor systems of interest. 
2.1 Action and field equations 

We consider a single light scalar field, and choose its action to be given by^ 



where a Weyl re-scaling is performed to go to the Einstein frame (which eliminates any (p- 
dependence from the Einstein-Hilbert term) and the scalar field is redefined to ensure its 
kinetic term is minimal (with a conventional factor of 2). Here G is the Einstein- frame 
gravitational constant, g^i,, is the Einstein- frame metric whose Ricci tensor is TZ^^, and S'm 
denotes the 'matter' action, involving all other observed fields (collectively denoted here by 



There are two physical choices made in writing this action, beyond the choice of using 
only a single scalar field. 

• We assume the absence of a scalar potential, which we assume is small enough to be 
negligible for the astrophysical applications of interest. This would necessarily be true 
if the scalar is relevant to cosmology, but is also the feature that quantum contribu- 
tions make most difficult to achieve in realistic models (unless there is an approximate 
symmetry, like shifts 4> ^ (p + c, for constant c). 

^Conventions: we use metric signature ( — h ++) and Weinberg's curvature conventions (differing from 
MTW only by an overall sign for the Riemann tensor). Units are chosen with h = c— 1. 
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We assume the matter action, Sm, only depends on (p and g^^i, through the one 'Jordan- 
frame' combination g^i^ = A?{(f))g^j^i,, where the conformal factor A{(j)) is a function 
whose form would be specified within any particular fundamental theory. This kind of 
coupling actually arises in specific models (such as if (j) arises as the breathing mode 
for the geometry of extra dimensions), and has the attractive feature that it naturally 
evades many of the strongest observational constraints on violations of the equivalence 
principle. 



The field equations obtained by varying (2T) are 

+ 2d^(pd^(l) + SvrG {t^^ - \Tg^,^ = (2.2) 

□(/> + 47rGa((/))r = , (2.3) 

where □ = S'^^'V^Vj/ is the d'Alembertian operator built using the Levi-Civita connection of 
g^y, while 

T/^- = --L^ (2.4) 
\r^og^y 

is the Einstein-frame energy-momentum tensor, T = g^'^T^y is its trace, and a((/)) = A'((/))/A((/)) 
defines the scalar-matter coupling function in terms of the function A{<j)). 

2.2 Observational constraints 

Because it is g^;^ that appears in the matter action, it is the geodesic of this Jordan-frame 
metric along which the trajectories of matter particles tend to move (in the absence of other 
forces). Upon taking the post-Newtonian limit of the Jordan- frame metric g^j^, one finds that 
the effective Jordan- frame gravitational constant, measured in asymptotic Einstein- frame 
units, ^ is 

G = G\l + a^{>Poo)\, (2.5) 

and that the Jordan-frame parameterized post-Newtonian quantities whose values differ from 
those of GR are 

2[1 + a^(0oo)J'' l + a^(0oo) 

where (j)oo is the asymptotic value of the scalar field far from the source and b{(p) = da(0)/d0 

i- 

Constraints on these PPN parameters from solar system observations provide the best 
bounds on the model, with I7 — 1| < 2.3 x 10~^ inferred from Cassini tracking, and |/3 — 1| < 

2.3 X 10~^ from lunar laser ranging This corresponds to the coupling bounds 

a\cPi')< 1.2. and a\<Pi^)\b{<Pi^)\ < 4.6 ■ W'^ (2.7) 

where denotes the asymptotic value of (j) as one leaves the solar system. 

^These are units in which the Einstein-frame metric is asymptotically Minkowski: diag(— 1, 1, 1, 1). If one 
instead uses units in which the Jordan-frame metric is asymptotically Minkowski, then G — Gyl^(<^oo)[l + 
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2.3 (Quasi) Brans/Dicke scalars 

The simplest scalar-tensor theory is Brans-Dicke theory |22|, for which a(0) = is constant. 



In this case - see eq. ( |2.7D - solar system tests constrain < 3.5 • 10^^, and so all of the 
predictions of Brans-Dicke theory are very close to those of GR. 

The next-simplest theory, which we call quasi-Brans/Dicke (qBD) theory, allows a{(j)) to 



vary slowly with (p [23|: 



A((/>) = exp(a,(/) + i , a{(l)) = as + bs<p . (2.8) 
This introduces a variety of new phenomena because it makes the strength of the scalar- 



matter couplings depend on (p, and so allows them to vary with position and time |24, 25]. 
This means that couplings in exotic environments (like stellar interiors) could be stronger than 
naively expected without running into conflict with the strong solar-system bounds mentioned 
above (this is similar in spirit to, but different in detail from, evading these bounds through 
matter-dependent scalar self-couplings ^ ^)- Assuming (j) is defined such that (?!)—)• 
asymptotically far from the Sun, eqs. (^) show that the strong bound on implies that 
solar system bounds do not strongly constrain bg. 

The constraints on bs are comparatively weak, and the best come from studies of binary 



pulsars [p^, 30, 31]. The precise timing measurements that are possible for binary pulsars allow 



their orbits to be accurately measured over long periods of time, and comparing measurements 



with the predictions of the qBD model leads to the constraint bs > —5 [10, 11, |l^, 13]. 
Measurements of the redshift of spectral lines from neutron stars leads to a weaker constraint 
65 > —9 [^]. The main uncertainties in these bounds come from the poor understanding of 
the nuclear equation of state appropriate for neutron star interiors. 

Observations disfavor negative bs because for bs < —4, compact objects like neutron stars 
exhibit a phenomenon called spontaneous scalarization JI8]. For sufficiently dense objects. 



whose precise threshold density depends on bs and the nuclear equation of state, it is possible 
for the star to support a nonzero scalar field even though the scalar coupling vanishes asymp- 
totically far from the star: a{(j)oo) = 0. Furthermore it is known that when scalarization takes 
place it is the stable solution to the field equations ]|3^, |l|. Scalarized neutron stars tend to 
be disfavored by binary pulsar observations because scalarization significantly changes the 
dynamics and radiation generated by a neutron star. 

Collapse processes involving scalarized neutron stars have been investigated by a number 
of authors J^], who found that the waveform of the emitted gravitational radiation depends 
strongly on bs- Moreover, scalar-tensor gravity allows monopole and dipole radiation, which 
are forbidden in GR. Thus, future measurements of gravitational waves may lead to improved 
constraints on 6, ]I]1 



3. Stellar structure 

One of the potential uncertainties when trying to constrain scalar-tensor models using astro- 
physical tests of gravity is the strength of the scalar field that is supported by objects like 
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the Sun or a neutron star. For macroscopic astrophysical objects made up of weakly coupled 
constituents one's intuition is that the scalar coupling to the macroscopic object should be 
proportional to the scalar coupling to each of its constituents, and we shall see in this section 
that this intuition is generally borne out for weakly coupled non-relativistic systems. We shall 
also see that it can fail for relativistic systems, even in the limit of weak scalar coupling. 

In order to do so, we next summarize how the scalar field alters the physics of stellar 
interiors, since it is the matching to this that dictates the properties of the external field 
configurations to which external observers have observational access. 



3.1 Equations of hydrostatic equilibrium 

The equations of hydrostatic equilibrium in scalar-tensor gravity were first studied in 



and 



I, n, n, 0, n, n, g, n, n, m, |7 



51], H, |5|, |5J 



subsequently in |3g, ^ 

Following [^] we model the stellar interior as a static, spherically-symmetric and perfect 
fluid in local thermal equilibrium, locally characterized (in the Jordan frame) by its pressure, 
P, and mass-energy density, p. Time-independence and spherical symmetry allow the use of 



Schwarzschild-like coordinates for the Einstein-frame metric interior to the star: 

dr^ 



gaS dx°dx^ 



+ 



+ r^dfi^ 



(3.1) 



1 - 2^(r) 

where dfJ^ = d9^ + sin^ 6 d(f>^ denotes the usual round angular metric on the 2-sphere and 
z^(r) and fj,{r) are to-be-determined functions that depend only on the radial coordinate r. 
The Jordan-frame energy-momentum tensor for matter is defined by 

2 6Sm 



and is related to the Einstein-frame energy-momentum tensor by 



nap 



T^p = A\<P)T^p, T^^ = A\(t>)Tj and T'^^ = A\<P)T 

where indices on T^p are raised and lowered using the Jordan- frame metric, Qaj^. 
Being a Jordan-frame perfect fluid, the energy-momentum tensor has the form 



(3.2) 



(3.3) 



(3.4) 



Tap = {P + P)UaUp + PgaP , 

where P = P{r) and p = p{r), and Ua is the Jordan- frame 4- velocity of the perfect fluid, 
given in co-moving coordinates by 

Ua = e''/^A{(t>)5i, (3.5) 

so that g°'^UaUf^ = —1. 

When writing the fleld equations it is convenient to scale out a dimensional factor of the 
central density, po := p(0), from the density and pressure, 

P(r) . , , p{r) 



p{r) :- 



and Q{r) 



po 



(3.6) 
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in terms of which the equation of state is specified by writing g^r) = Q\p{r);po]- Here 

^(0) Po 



Po 



(3.7) 



labels the star's central pressure in the same units, and the po-dependence of ^>[p;po] is meant 
to emphasize that the functional form of the equation of state, P{p), in general changes as 
one changes the central density. Notice that our notation implies the identity q{po',Po) = 1- 

With these choices the Einstein field equation, eq. ( |2.2| ), in this geometry boils down to 
the following three conditions: 

^2 



rf,' + fi = A7rGpor^A\<P) q{p) + y (1 - 2^)0" ■ 



P 



-[q{p)+p] 



A-KGpor^A*{4>)p + 
r(l - 2/i) 



+ 



+ a 



87TGpor'A^{^)p + 2p , 2 
; ^ \- ' 



(3.8) 
(3.9) 

(3.10) 



where primes denote derivatives with respect to r. The scalar wave equation, eq. (2.3), 
similarly becomes 



1 - 2p 



a{(l))[0{p) - 3p] + r(p'[g{p) - p] 



r{l-2pf 



(3.11) 



These equations are to be integrated subject to the following initial conditions at the centre 
of the star: 

p{0) = 0, p{0)=po, m = <Po, </.'(0) = 0. 
Writing eqs. and (|3lO| ) as 



2p' + {p+Q) u' + 2a{^)(p' =0, 
and integrating once allows v to be written in terms of p and (p: 

V = -2f{p;po) -2ln A{(j)) + const , 

where 



(3.12) 
(3.13) 



fip;po) 



dp 



(3.14) 



(3.15) 



/po P + e{p) 

Regarding eq. ( |3.8D as a linear, first-order differential equation for the product rp, allows 
it to be solved to obtain 

l4/ 



1 



P 







+ 4t^GpoA\4>) Q{p) 



(3.16) 



where the integration constant is chosen such that ( p.l2| ) holds. This expression shows that p 
is always non-negative. It also shows that p{r) should not be interpreted as the mass-energy 
inside the ball of radius r, unlike in GR. 
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In principle, eqs. ( [3.8D through ( 3.11 ) can be integrated numerically, starting at r = 
and continuing out to larger r. In practice, this system of equations is singular at r = 0, so 
numerical integration must be started at some small positive r = rg. The initial conditions 
at ro can be obtained from the power series expansions that are dictated by the equations of 
motion themselves: 

^J'{r) = ^ -r +0{r ), 

p{r) = Po + ^^^3°^° (P0 + 1) [ali^Po - 1) - (3po + l)V + 0{r^) , 



where = ^('?^o) ^"^^ ^^o = o(i?i'o) 



ao(3po - \y + O(r^) , (3.17) 



3.2 Matching to exterior observables 

The integration within the stellar interior continues until eventually the pressure p vanishes. 
The value r = R where this happens defines the (Schwarzschild coordinate) radius of the star, 
beyond which the appropriate solution instead solves the 'matter- vacuum' field equations with 
p = P = 0. Of course, for generic scalar-tensor theories it might happen that p never actually 
vanishes, since unlike for GR p need not be a monotonically decreasing function. We comment 
where we can on the stability of these configurations below, although as we shall see they are 
unlikely to happen sufficiently close to the Brans-Dicke limit and for non-relativistic equations 
of state. 



Exterior solutions 

For r > i?, an exterior solution is required to satisfy the vacuum field equations. 



and = , 



(3.18) 



for which a closed-form static and spherically symmetric solution may be found by taking 
(f) = (j){r) and using the metric 



-e^" df + e-2« dr2 + e^'" 



d6'2 + sin^ e dip' 



(3.19) 



where u = u{r) and v = v{r). The field equations for (p, u and v have the following solutions 



1 - - 

r 



1-x / 

rUl--j and e^* = e^^^- ( ^ " " ) ' (^-^O) 



where the real constants x and q satisfy x'^ + = 1 and are otherwise arbitrary, leaving three 
free integration constants: x, i and (poo- 
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These three constants are more conveniently rewritten in terms of the large-r hmit of the 
solution when expressed in the original Schwarzschild coordinates, 



.21 2 GM , , 2 

ds2 = -1 + + ... df + 



1-^ + ...), (3.21) 



where M is the system's ADM mass |^6[| in the Einstein frame and Q defines its 'scalar 
charge'. For systems where the scalar coupling, a{(p), is field-dependent, the asymptotic value 
of the scalar field, (p^, may usefully be traded for the asymptotic value of the scalar coupling 
strength, Ooo = a((/)oo)- 

Together with the star's radius, we see that from the point of view of an external observer 
there are four independent bulk parameters that characterize any such a star in scalar-tensor 
theory: M, Q, R and (/jqo- These can be calculated in principle as functions of the stellar 
equation of state by matching the exterior solution to that of the interior at r = i?, implying 
that they are functions of the two parameters, (po and po, that define the initial conditions 
of integration at the stellar center, eq. ( 3.12| ). One of our main goals is to identify the two 



relations that must hold amongst the four external parameters (for any given equation of 
state), 

Ci{R, M, Q, </<oo) = UR, M, Q, </<oo) = , (3.22) 

generalizing the familiar mass-radius relation, R = R{M), that expresses the content of 
stellar structure within GR. The explicit form of these constraints for scalar-tensor models is 



discussed in section 4.2 below. 

Physically, we expect the value of (j)oo to depend on physics external to the star or stellar 
system of interest, governed by the local properties of the galaxy in which the stars are 
located, (poo could also depend on cosmological time if (j) is light enough to be evolving over 
cosmological time intervals. Thus, only one combination of the four parameters M, Q, R and 
(j)oo can vary in practice from star to star within a specific galactic neighborhood at a given 
cosmological epoch. In particular, this means that the scalar charge Q and mass M of a star 
are not independent parameters for any stars within a local neighborhood at a given epoch. 

This result greatly simplifies the phenomenological analysis of binary pulsars. Instead 
of having two independent parameters describing each of the two stars in the binary system, 
there is in practice only one and so it suffices to describe the observational constraints as a 
function of the masses of the two stars, just as in GR. In principle the dependence on (p^o 
could complicate the combining of the implications of many pulsars that are located far from 
one another, but a simple estimate shows that (j) does not vary strongly across the galaxy, so 
in practice (j)oo can be taken to have a common for all stars in the galaxy. This approximation 
is implicitly used in refs. fl^, 0, |l2|, when combined data from several binary pulsars are 
plotted on one theory-space exclusion plot in terms only of the masses of the two components. 
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To estimate the variance of (p{r) across the galaxy, we may use the non-relativistic hmit 
for which g'^'^T^^ = A'^{(p)g''l^fo,f} = ^^(0)(3p - p) ~ -A^icp) p and so 



ATTGa{4>)A^{(t))p. 



(3.23) 



If, on the right-hand-side, we assume < 



to be approximately constant, then a{(j))A^{(l)) 



cigA'i is also constant. Then the only r-dependence on this side comes from p, and we may 



estimate how large a deviation from constant (j) is implied by eq. ( 3.23 ). For r large enough 
that p is dominated by Dark Matter, the density profile is p ~ (i/r)'^ so that orbital 
velocities are r-independent: Vg ~ GM{r)/r ~ AnGpgi'^. Then eq. ( 3.23| ) implies 

r 



(r') + inGagA^nglHn 



anAtvl In 



(3.24) 



This is clearly small because it depends only logarithmically on r, and both ag and Vg are 
small. See Appendix B of p^] for a related discussion. 



Matching at r = R 

Requiring continuity of the exterior and interior profiles across r 
external parameters M, Q, and (pao must satisfy 111 



GM 



IC 



s :- 



exp 



R 2V1 - 2// 
Q_ _ 2R(j)'^{l - 2p 
M ~ 



IC 



— arctanh ( — 



C 



C 



+ 



JC 

2R(l)'^{l-2p, 
C 



arctanh 



J 



J 



R implies the three 

(3.25) 

(3.26) 
(3.27) 



where 



J 
K 
C 



2(l-/i^) + i?Vv?(l-2pO, 
2p^ + R^(l)'^{l-2p^), 

V4/i2 + AR^'^il - p,){l - 2p,) + i?Vl^(l - 2p^Y 



(3.28) 
(3.29) 
(3.30) 



and p-i, := p{R), 0^ := 4'{R), and cp'^ := (p'iR). 

It is sometimes useful to decompose the ADM mass as M = Mg + AM, where AM is 
the gravitational binding energy and Mg is the baryonic mass, defined by 



^ , nA^(6)r'^ 
dr 



Mg := inmo / dr ^iL. , (3.31) 
Jo VI - 

where n(r) is the local baryon-number density, and rrih is the average mass of an individual 
baryon. This baryonic mass can be related to the stellar pressure and density if the star has 
constant entropy per nucleon, since in this case it follows from energy conservation that 



n 



n 



0. 



(3.32) 
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This equation can be used to write n in terms of p, 



n 



no 



P + Qjp) 



3-/(p;po) 



(3.33) 



where no = n(0) and f{p',Po) is defined in equation ( p.l5| ). The baryonic mass then becomes 



R 



p + q{p) 



A^{(j))r'^dr 



(3.34) 



The quantities s and a a given by the matching equations, ( 3.25 ) and ( p. 26 ), are in 
themselves useful because each has a physical interpretation, s is called the self-gravity^, 
or compactness, of the star and it provides a dimensionless measure of how relativistic its 
gravitational field is at the stellar surface, r = R. For non-relativistic stars s <C 1, and 
in general relativity s < | for any star within which the mass-energy density, p, is a non- 
increasing function of r, a result known as Buchdahl's theorem [f!g|].^ 

Similarly, the quantity Qa can often be interpreted as the effective scalar-matter coupling 
of the star as seen by the observer at infinity. This can be understood by considering the 
lowest-order non-relativistic interaction energy between two widely separated stars, A and B, 



GMaMb GQaQe 



G.nM.M^ 



rAi 



(3.35) 



where 

Gab = G{1 + ttAag) (3.36) 

is the effective Jordan-frame gravitational constant between the two stars, and r^^ is the 
distance between them. All quantities are measured in units in which the Einstein-frame 
metric is asymptotically Minkowski^. 

The connection between and the scalar coupling can also be seen in another way. 
Since in the non-relativistic limit we have M ~ Mg, and since each individual baryon if 
separated to infinity would microscopically satisfy qb/rrib = aoo = a(0oo); we expect that for 
macroscopic non-relativistic systems Q/M = {Nqt,)/{Nmb) = Ooo as well, so 



lim Ua = a( 
s-s>0 



and 



lim Gab 

S-5>0 



G, 



(3.37) 



where G is the Jordan- frame gravitational constant defined in equation ( p.5| ) [Q. We shall see 
below that the prediction Q/M = a{(poo) is given by a more detailed matching calculation, 
where it can be seen to hold independent of the stellar equation of state in the limit of both 
non-relativistic sources and weak scalar coupling. 



^Note that Will |^] uses s to denote the sensitivity of the mass to the gravitational constant, which he 
defines as — (91og Ai")/(9 log G). This definition of s differs from our definition, but has the same order of 
magnitude. 

''The potential generalization of this result to scalar-tensor gravity is investigated in js^]. 

^As mentioned earlier, to convert to units with an asymptotically Minkowski Jordan-frame metric multiply 



Gab by A' 



multiply Tab by A{<t>oo), and divide Ma, Mb and Uab by A{(l>a 
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The black hole limit 



For sufficiently massive stars the only stable configuration in GR is a black hole, for which 
there is only the exterior solution, eq. ( 3.2[1| ) with x = 1 and q = [^, 21 1. In this case the role 
of matching to an interior geometry across r = R becomes replaced by the boundary condition 
that the geometry not have a physical singularity at the event horizon, r = Rbh, making 
this the natural equivalent of the stellar radius, R, for a black hole. Since for spherically 
symmetric stars in GR we have R^h = 2 GM, we see that for black holes the mass-radius 
condition predicted is particularly simple: s = GM/Rbh = | > |- 

For a black hole coupled to scalar fields the exterior solution, ( 3.20|) , always has a curva- 
ture singularity at r = £ whenever q This shows that only q = describes a legitimate 
black hole, for which (p = cp^ is constant and arbitrary (this is unchanged if Ooo = o-ifpoo) 0) 
because there is no matter exterior to the black hole with which to couple). The metric is 
then given by the Schwarzschild geometry, just as for GR, implying s = ^ and Q = (for all 
0oo) even when scalars are present. This is a special case of the no-hair theorems, originally 
formulated for gravity in |59] and extended to scalar fields in [BO]. 



Stability 

Once a solution of equations (^)-( 3.11 ) 
is obtained, and the external parameters 
( |3.25| )-( |337| ) are calculated, it must be checked 
whether the solution is stable against per- 
turbations. To do this properly requires 
working with the time-dependent equations 
of stellar structure. We here instead follow 
[ p^ , |3^ , |l] and perform a simplified analy- 
sis. 

The main idea builds on the observa- 
tion that (for general relativity) the equa- 
tions of hydrostatic equilibrium are equiv- 
alent to the problem of extremizing M, 




Figure 1: Mb versus M for constant-density stars 
in Brans-Dicke theory, with 0oo — and aj. = 0.2. 



while keeping Mg fixed 1 2C ] . The idea used 

in 1 18, |3^, P is to assume the same is true for scalar-tensor theories (with (p^o fixed as well as 
Mg) although we have been unable to prove this so far in scalar-tensor gravity. 

For given values of Mg and (poo , one generically finds that there can exist multiple stellar 
configurations with different values of M. This is illustrated by the plot of Mg versus M 
given in figure I for a one-dimensional family of stellar configurations that share a common 
value for <poo- Whenever two or more values of M are are obtained for a given (poo we take 
the configuration with the lowest value of M to be stable, and the others to be unstable. 

3.3 Stellar structure for quasi-Brans/Dicke scalars 

Of course all of the predictions for quantities like Q/M depend in principle on the details of 
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the coupling function, a(0), in addition to depending on the stellar equation of state. In this 
section we specialize the equations of the previous section to the quasi-Brans/Dicke (qBD) 
model of eq. (|2.8|) , since this boils the dependence on a{(f>) = as + bs4> down to the dependence 
on the two parameters and bs- 

Because the gravity-scalar part of the theory is invariant under constant shifts, </>—>■ 0+c, 
this transformation can be used to set to zero with no loss of generality, provided that 
bg 7^ 0. With this choice, the matter coupling can be seen to have a reflection symmetry 
(j) —7- —(j). We henceforth choose this convention for a^, making it sufficient to follow the 
dependence of observables on 6^- 

Consider now integrating the field equations starting from r = 0. The information that 
used to reside in now resides in the value of the field (or coupling) at the stellar origin, 
ao := a(^o) = ^s'^o- Suppose first that the initial value, cpo, chosen for (j) at the stellar center is 
such that the coupling vanishes there, a{(j3o) = 0. In this case the reflection symmetry implies 
(p{r) = for all r is a solution to the field equations, and so a{(j)) = oq = 0, everywhere within 
the stellar interior. (That this is a solution can be seen by direct inspection of eq. ( ^.11 ).) In 
this case the equations of stellar structure, ( |3.8D through ( p. 11 ), uniquely reduce to those of 
GR. 

If instead oq = a(0o) = &s0o 7^ 0, then, after the change of variables 



4 „2 



(/9 = (0 - 0o)/ao , u = SttGpqAq r 



the equations of stellar structure simplify to 



P 



2u 



■){p) e4«o</'(i+''-'^/2) 



-[q{p) +p] 

(3 - Ap)^ 



2m(1 - 2//) 
+ 



+ 



f 00^(1 — 2^)99 



•2 



+ 0992(1 + ULP + bgip) 



g4ag</j(l+6s(/3/2) 



4(1 - 2fi) 
1 + bsip)[g{p) - 3p] + 2[q{p) -p]uip 



(3.38) 

(3.39) 
(3.40) 

(3.41) 



2n(l-2/i) 8m(1-2^) 

where dots denote derivatives with respect to u. These equations show that ao = a(</'o) = 
a{r = 0) only enters physical observables through its square: Oq. In terms of these variables 
the initial conditions are 

MO) = 0, 



p(0)=po, <^(0)=0, (^(0) 



1 

12 



(1 - 3po) . 



(3.42) 



The profiles /i(u), p{u) and ip{u) obtained by integrating these equations depend only on the 
three parameters Oq, bg, pQ, as well as on the choice of equation of state g{p). 
The baryonic mass, eq. ( 3.34| ), then becomes 



27rmono 



27rmono 



u 



P+Qjp) f(p;po) 3a^Ml+bsV/2) 



(8^Gpo)3/2Ai] 



dn — e 

M{al,bs,po) , 



1 - 2fi 



(3.43) 
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where 

U{alh,p^) = ST:GpoAlR^ (3.44) 

is the value U = u{R) corresponding to the boundary of the star. For a given equation of 
state, U depends only on Oq, 6s and po, because it can be computed by finding the point at 
which p{u) vanishes. 

Note that the quantities J , K., and C, defined in equations ( 3.2^ ) through ( 3.30| ), depend 



only on fj.^ = fi{U) and R'^(j)''^{R) = Aa^U'^tp'^iU). Therefore, the compactness s — given by 
equation ( |3.25| ) — also depends only on Oq, 6s and pq. The matching conditions at the stellar 
surface, eqs. (|33^ ) - ( |337| ), can be re-written using these variables as 



^{U) + 2^--) arctanh ] := ^(ag, 6s,po) , (3.45) 



ao C \J 

and 

= = ^ := aoA{aQ, bs,po) , (3.46) 

which define the quantities and A. 

Properties of the pressure profile 

This section briefiy pauses to investigate whether the pressure, p, decreases monotonically, 
and whether the equation p{u) = must have a solution. 

It is clear that (unlike for GR) dP/dr can be positive for sufficiently large couplings and 
pressures. This follows from evaluating equations ( |3.40| ) and ( ^.42 ) at the stellar center, 



P(0) = -^^[l + 3po + ag(l-3po)] , (3.47) 
which (for pQ > 0) is strictly non-positive unless Po = Po/po > I and 

2^2 3po + 1 fn Aa\ 

Oo > «Ocrit ■= 3p^_l ■ (^•^^) 

In general, to determine whether the pressure profile is monotonically decreasing everywhere, 
it is not sufficient to look only at p{0). In the qBD theory with 6s > 0, it can happen that 
p{0) > 0, p{Ui,) = for some > 0, and p{U) = for some U > u^. Such a star has its 
maximum pressure in between the centre and the surface. Similarly, if 6s < 0, it can happen 
that p(0) < 0, p{u^) = for some > 0, p(u) > for all n > 0, and p{u) — )• oo as u — )• oo. 
Such a solution describes an object of infinite extent, and seems unphysical. 

By contrast, for Brans Dicke theory our numerical calculations suggest that this kind of 
complicated behaviour does not take place, and the pressure profile p{u) cannot have any 
local extrema. We have been unable to prove this analytically, but we have shown that if 
Og = Oocrit' then Il{u) = t] for all u, for all equations of state. This constant-pressure solution 
appears to be the boundary between the solutions with p > 0, and those with p < 0. For 
Brans Dicke theory Salmona |^4| has shown that if p/g < 1/3 everywhere, then p decreases 
monotonically everywhere. 
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If we assume that this simple behaviom' of the pressure profile in Brans-Dicke theory is 
correct, then imposing certain requirements on the properties of the solutions of the equations 
of hydrostatic equilibrium leads to constraints on Uq. If we require that the pressure is 
decreasing for all pQ, then it follows that Oq < 1. 

A more conservative constraint on Oq can be obtained by using some information about 
the equation of state. For these purposes it suffices to consider only neutron stars, because 
only these are relativistic enough to have po > 1/3. But neutron-star interiors can be modelled 
as relativistic polytropes with 

p Ku /n\i , P /n\"/ , , , 

— =n+ - and — = /tn - , (3.49) 

nib 7 ~ 1 i^b ^nJ 

where n is the baryon number density, n = 0.1 fm~^ is a typical nuclear density, rrn, = 
1.66 • 10-2^ g is the mass of an average baryon, and /t and 7 are constants (7 is called the 
polytropic index) . Notice that these choices imply the central density and pressure are related 
by po/rub = no + Kh{no/hy /{'y — 1), and so 

1 = ^ + ^, (3.50) 

Po 7-1 

implying the maximum value of po that is possible is p^^^ = 7 — 1. The functions q{p) and 
f{p) are similarly given by 



Qip;po) 



(3.52) 



fip;po) = In 



7-1 

7-1 -po + jpipo/p)'^^'^ 



(po + l)(7- 1) 



Requiring that the pressure is decreasing for all relativistic polytropes yields the con- 
straint 

<.§<|^. (3.53) 

For the neutron-star equations of state EOS II and EOS A of ref. [^, these bounds evaluate 
to |ao| < 1.29 and |ao| < 1.26, respectively. Although nowhere near as good as the solar 
system bounds of eqs. (^), they are complementary because they apply to the coupling deep 
within a neutron star interior, and rely only on general considerations, rather than precise 
observations. 



3.4 Non-relativistic limit 

In most stars, relativistic effects are negligible, allowing us to take p = P/po <^ 1 and /i ^ 1. 
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In this case, the equations of stehar structure, eqs. (|3.8|) through ( 3.11 ), simphfy to 



r/i' + /i ~ ATTGr'^A'^{4))p + 



P' 



fi r _ 



12 



+ a 



(f)" ~ AttGA^{(I))p a(0) + rcj)' 



24)' 
r 



r 



(3.54) 
(3.55) 

(3.56) 
(3.57) 



where the energy density, p, is equivalent to the mass density in the non-relativsitic hmit. 
The matching conditions, eqs. ( 3.25| ) through ( |3.27 ), similarly simplify to 



GM 



p{R) + 



(3.58) 

(3.59) 
(3.60) 



GQ 
R ' 



(3.61) 



0OO = (l){R) + R<p'{R). 
Notice that these last two matching conditions quite generally imply 

<t){R) = (t>o 

as is appropriate for the non-relativistic limit of the known external solutions. 
Newtonian polytropes 

The equation of state that is of most interest for Newtonian systems is that of a polytrope, 

P = Kp^+^'^ , (3.62) 

where x is the polytropic index (a constant that need not be an integer), and X is a constant. 
We briefly specialize to this equation of state here for later use in subsequent sections. 

Specializing eqs. (|3.54|) - (|3.57|) to this equation of state, ( 3.62| ), the above equations 
simplify after changing to dimensionless variables 



rsW , 



(3.63) 



where po = p(0), and is a length scale that is to be specified later. Then equations ( p. 541 ) 
- (@) b 



become 



dw 



-— + GA^{(t))e^ \a{<i))+w4) 
w L 

-GA\<P)wH^-^4>\ 



-w 



+ 17 + a{<t))^ 



(3.64) 
(3.65) 
(3.66) 
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where dots now denote d/dw, and the dimensionful parameters are all rolled into the new 
dimensionless constants, 

C := K{x + l)Po^^ , C := A^Grlp^ . (3.67) 
In these variables the initial conditions are 



0(0) = 1, ^(0) = 0, 
and (/)(0) = , 0(0) = . 

We now specialize to the quasi-Brans/Dicke models of eq. ( |2.8D , for which we had 



al = a^((^o) = [a-s + f's'/'o)^ and ip 



4>- (Po 



ao 



If we choose 



so that 



Al y 4ttGpo Al V 4ttG 



1 ./ -^(X + I) j-X)/2X 



then equations ( f3.64 ) - ( 3.66| ) become 



C = CA 



) 



dw 

and p = -w (^6 + ^agW*^^ + ao(l + &s'/')(^^ , 



and the initial conditions of eq. ( p.6§| ) for (p now are: 

ip{<d) = , ip{<d) = . (3.76) 
This system implies the solutions have the following power-series expansions near w = 

(3.77) 



e{w) = l-^{l + al)w^ + 0{w^) 
6 



When written in terms of the variables 9, and w, the matching equations, (^ 
( |3.6C ), become 



wUe{W) +al[l + hs^{W)]ip{W) 



GM 

and (t)^ = (t)Q + ao[^{W) + Wp{W)] , 
where W = w(i?) = R/rg denotes the stellar boundary. 



;3.68) 
3.69) 

3.70) 

3.71) 
3.72) 

3.73) 
3.74) 
3.75) 



3.78) 



53) - 



3.79) 

3.80) 
3.81) 
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4. Solutions for weak central coupling 

Most of what is known about the solutions to the equiUbrium equations derived above is based 
on integrating them numericahy, reveahng several surprising features such as the phenomenon 
of spontaneous scalarization [Q. But the regime of most practical interest is weak coupling, 
Cq ^ 1, and since interesting phenomena like scalarization are already present in this limit, it 
is worth exploiting the simplicity of the weak-coupling regime at the outset, both to simplify 
the numerics required and (in some cases, see below) to allow analytic solutions to be obtained. 

Our goal in this section is to systematically expand in powers of the scalar coupling at the 
stellar centre, Og = a^{4>o) = a^(r = 0) ^ 1. By comparing these perturbative results with 
direct numerical integrations, we show that for small Oq their domain of validity typically 
covers the entire stellar interior. 

Our motivation for pursuing the simplifications introduced by this expansion is the ease of 
generalizing to different kinds of scalar couplings and to different equations of state. However 
in this paper we confine our attention to the well-studied qBD case, a{(p) = + bgCp, in order 
to better compare with known results. 



4.1 The weak-central-coupling expansion 



The weak central-coupling expansion is clearest in the case of qBD models, for which the entire 
coupling function, a{(j)), is determined by the two parameters ao = a(</'o) and bg. In this case 
the weak central-coupling solutions are obtained by expanding the differential equations in 
powers of Og. For Brans-Dicke theory (bg = 0) the coupling is constant and known to be 
small, al = a'i< 1.2 x 10^^ 

For qBD theories the lowest-order expansion of the equilibrium equations, ( 3.39| ) through 
( |3^ ) gives 



p + 



2u 4 
[q{p) +p](2/i + -up) 
4«(1 - 2fj,) 

(3 - 4^)vi 



0{p)^ ( 1 + ^ ) + n(l - 2^)vi2 



- -ao[Qip)+P] 
1 



+ ip{l + u(p + bsf) 



2m(1 - 2fi) 8m(1 - 2fi) L 

_ ag(/?(l + bsip/2) 
2u{l - 2/x) 



1-2// 

{I + bsip)[g{p) -3p]+2[g{p) - p]uip 



(1 + bs^Mp) - 3p] + 2[0{p) -p]u^ 



(4.1) 

+ 0{4) (4.2) 
(4.3) 
+ 0(a' 



A), 



where dots denote derivatives with respect to u = SttCpqAqt'^. Notice that the leading 
contribution to the 93 equation depends on the self-coupling bg, even if Cg — )• 0. The boundary 
conditions are as before: /i(0) = 999 = 0, p{0) = po and (p{0) = j2 ~ ^Po) (and so ip{0) > 
provided po < | ) . 
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We seek interior profiles p{u) and (p{u) obtained by integrating these equations 

subject to the series ansdtze, 



= ^ (n) al' , p{u) = ^ {u) ag* and ip{u) = ^ V3(i) {u) a^' , (4.4) 

i=0 i=0 i=0 

with the leading expressions for pQ (u) and no (u) agreeing with the results from GR. In partic- 
ular, because Og is small, the pressure profile decreases monotonically, ensuring the existence 
of a solution for R of p{R) = 0. Because of the explicit factor of oq appearing in the definition 
if := ((/> — (/>o)/ao, given a solution for ^, p and ip correct to order a^^, we have a solution for 
(f) that is valid to order Oq'^'*'^. Such solutions are obtained explicity for k = and k = 1 for 



incompressible stars in sections 5.2 and 5.3 below 



4.2 Perturbative relations amongst observables 

This same Oq expansion is inherited by the expressions relating the external physical parame- 
ters, M, R, Q and (poo, by virtue of the matching conditions at r = i2. Although it is difficult 
to characterize these constraints analytically in the general case, an expansion in powers of 
Oq = a^(</>o) allows some progress to be made. To this end write 

U{ao,bs,po) ^ C/(o)(po) + ^(i)(&s,?'o) 00 + 0(00) 
J"(ao,6s,Po) - J^{o){bs,Po) + J^{i){bs,po) al + 0{aQ) 
A{ao,bs,Po) ^ A^o){bs,Po) + A(^i){bs,po) al + 0{a'^) (4.5) 
s{ao,bs,po) ^ S(o)(po) + S{i)(&s,Po) 00 + 0(00) 
M{ao,bs,po) ^ M(^o)iPo) + M(i){bs,po) al + 0{4) . 

where s = GM/R, T = {(poo — (po)/o-o, A = Q/Mao = Oa/oo and so on. 

We seek two constraints among the four quantities M, Q, R and (poo, ^-nd it is convenient 
to write the first of these as a relationship between M, R and c; 



, and the second 

relationship between = Q/M, s = GM/R and Ooo = a{(poo)- The convenience of this choice 
comes from the GR limit, for which the first constraint becomes the usual M-R relation, and 
the second constraint degenerates into something vacuous: = 0. 

An important point about these constraints is that the dependence of observables like 
M, Q and R (or U) ) — arises completely through their dependence 

on ao = cl{(Po)- Obtaining this dependence therefore requires a relation between the scalar 
field at the origin and infinity: (po{(poo)- This is accomplished by using the function J^, whose 
definition - see eq. ( 3.45 ) ~ states (poo = </'o + oq J^, and so 

A{(Poo) = A{(Po) exp \al (l + ^] , (4.6) 

(4.7) 



a{(Po){l + b,T) 



It is tempting to ask at this point whether we are working too hard. In particular, since 
it is Ooo and not ao that directly controls the strengths of interactions that we see, being 
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asymptotic observers, perhaps we could avoid the exercise of trading oq for a^o by directly 
expanding the field equations in powers of a^o rather than oq- The reason we do not do 
so — and indeed the point of expanding in powers of ag — is that the mapping defined by 



eq. (4.7) between oq and Ooo is in general not one-to-one. This is the lesson of scalarization, 
which relies on Ooo = corresponding to several choices: ao = and ao ^ 0. It is the option 
of having a second choice that allows the star to support a scalar field {Q 7^ 0) despite the 
vanishing of Coo • It is the fact that integration of the field equations makes stellar properties 
single-valued in oq that makes this the natural expansion parameter. The existence of several 
branches to the function ao(aoo) means that stellar properties need not also be analytic in 
floo- We describe the relevance of this to scalarization in more detail below. 

The generalized mass-radius relation 

To obtain the leading form of the constraint generalizing the M{R) relation of GR, set Og = 



in equations ( 3. 391) - (|3.4[l|) . The result implies that ;U(o) and p(o) "^o '^ot depend on bg 



Consequently C/(o)(po)) which is defined by P(o)(f^(o))Po) = 0, also cannot depend on 6^ — a 
fact already indicated in eqs. (|4.5|). The same is then true for the compactness, 



s{U,cPoo;bs) = ^^s^o)[U^o),po{U^o))] + 0{al), (4.8) 

implying this constraint goes over to the GR limit to leading order in Og, even though the 
profile for ip{u) need not be trivial (for nonzero bg). Thus is reproduced the usual M vs R 
(or U) relation of GR. 

Figure m illustrates how the generalized mass-radius relation depends on in Brans-Dicke 
theory {i.e. bg = 0) for incompressible stars (discussed in more detail in the next section) and 
relativistic polytrope models of neutron stars. Each curve traces the relationship between 
M and R as po is varied, beginning at M = i? = where po — t- 0, and terminating at 
the point where po — >• p^^^. For incompressible stars, p™^^ is infinite. In general relativity, 
M (X R^, and the maximum M and R that can be supported against gravitational collapse 
are attained in the ultra-relativistic limit. However, once the scalar-matter coupling is turned 
on, the maximum values of M and R are attained at a finite value of pQ. As increases, the 
maximum values of M and R decrease. 

The equations of state EOS II and EOS A are defined in reference They are relativis- 
tic polytropes, with a maximum central pressure of p^^^ = 7 — 1, where 7 is the polytropic 
index. Their M — R curves are more complicated than those of incompressible stars, and 
turning on a weak scalar-matter coupling slightly shifts these curves. 



Numerically carrying out the stability analysis described in section 3.2 shows that in all 
cases, the stellar configurations become unstable after the first turning point where dM/dR = 
0. Thus, the scalar field destabilizes ultra-relativistic incompressible stars. 



For non-relativistic systems the matching condition simplifies to eq. (3.58), 



'(0) = M{o)(^(o);Po) < 1, (4.9) 
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and although numerical methods are usually required to follow the dependence on po{R), 
more explicit statements about the scalar corrections to this relation are possible for specific 
choices of equation of state. The examples of Newtonian polytropes and incompressible stars 
are considered more explicitly below. 

Scalar-coupling constraint 

For a given equation of state the second observable constraint, eq. ( |3.46| ), gives an expression 
for a_4 = Q/M in terms of oq, bg and pq. Equation ( |4.7| ) can be used to relate a(i;^oo) with 
a{(po), and on expansion yields 

a(0oo) ^ (1 + bsT(^o))ao + &s-^(i)«o + O(«o) • (4-10) 

Therefore, 

= § = aoAialM,po) ^ ^^^^ + 0{al) . (4.11) 

This expression diverges when 1 + 6s-^(o) ~^ 0- this limit, one must include the ©(oq) 
terms in order to obtain a meaningful result. This will be described in the section below. 

The dependence of a^/ao on s is shown in Figure ^ for several choices of oq = in 
the special case of Brans-Dicke theory (bg = 0), using for illustration an incompressible star 
equation of state (see §5, below). Each curve can be parametrized by pq. The starting point 
is at s = 0, ^ = 1 when pQ — )• 0, and the endpoint is reached in the limit pq — )• p^^^. 
These curves show that the small-ao limit works well in this case even out to very relativistic 
stars. Notice that a a is in general smaller than oq, with the suppression increasing for more 
relativistic stars. 

Fig. dH) also compares the amount of this suppression for several other choices for the 
equation of state, indicating that the suppression for relativistic polytropes (more about 
which later) is under-estimated for incompressible stars, although they all agree in the non- 
relativistic limit (for which S(o) ^ 1)) agreement with intuition. 

For incompressible stars a a eventually becomes negative. As seen in the insert in Fig. (^, 
this change in sign takes place before the onset of instability, so there exist stable incompress- 
ible stars with < 0. The scalar interaction between two stars A and B is attractive if 
a^Os > 0, and repulsive if CaCLb < 0. If at least one of A, B is an incompressible star, then 
both cases are possible. However, negative values of Oa are not seen for the more realistic 
equations of state, although these do approach = when extremely relativistic. This 
should be compared with the corresponding universal result, (5 = 0, found above for a static, 
spherically symmetric black hole. 

Spontaneous Scalarization 

In the above section, the scalar-coupling constraint was expanded in powers of oq, and it was 
found that at leading order, a a is a single- valued function of Ooo- In this section, we take the 
expansion to next order, and demonstrate that a a becomes a multi- valued function Ooo- As 
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a consequence, oq and can both be nonzero even when Ooo vanishes. The phenomenon 
where Q/M is nonzero even though a^o = is cahed spontaneous scalarization. To simphfy 
notation, write 

aoc = diaQ + d2al + 0{al), (4.12) 

aA = eiao + e2al + 0{al). (4.13) 

Dropping terms of order Og and inverting equation ( |4.12| ) yields 

ao = ujC+ + CjC- , (4.14) 

where 

C,= ,fe±V5. «=feV+ (4.15) 



2d2 \2d2J V3d2, 

and CO = 1, -e~*''/^ -e'""/^. Thus, 

= {ujC+ + CdC-)ei + (u;C+ + idC-fe2 + ©(a^) . (4.16) 

If > 0, then the a; = 1 solution is real, and the other two solutions are complex. Thus, 
specification of Ooo determines a unique stellar configuration. 

If D < 0, then all solutions are real. Thus, specification of Cqo determines three stellar 
configurations. 

In the limit Cqo 0, the u = 1 solution vanishes, and the other two solutions become 



<.„ = ±f|. (4.17) 

Therefore, whenever the quantity inside the square root is positive, there exist stellar config- 
urations with ttoo = and 




^Y\A+0{al). (4.18) 



d2j 

This is precisely the phenomenon of spontaneous scalarization. 

The coefficients di and are functions of bs and po, and they also depend on the equation 
of state. The solution obtained in section 5^ for incompressible stars can be used to calculate 

di = HeunG(a,g;a,/3,7, J;Z) - ^ "^J^" log ( } "^f" ) HeunG'(a, g; a, /3, 7, J; Z) ,(4.19) 

l + 3po \l + 3poJ 

^1 = A , HeunG'(a, q; a, ^,^,6;Z), (4.20) 

5^(1 + 3po) 

where the arguments inside the Heun functions are given in equations (5.34)-( |5.35| ), and 
Z = po/{3pQ + 1). Note that in the limit po — )• 0, we have di — J- 1 and ei — )• 1. In principle, 
the next coefficients d2 and 62 can be calculated using the solution found in section |5.3| . 
However, the resulting expressions are very complicated. 
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It follows from equation ( 4.17| ) that the onset of scalarization implies di = 0. Thus, 
all the points {bs,Po) at which scalarization starts can be found by doing a root search of 
equation ( 4.19 ). 

In figure (^), equation ( 4.19 ) is plotted versus po, for various negative values of bg. In 
all cases, di is a convex function of po, with one global minimum. If bg G (—4.329,0), then 
di never crosses zero, and there is no scalarization. If bg € (—oo, —4.329), then di has two 
crossings of zero, which correspond to the scalarization which has been extensively studied in 
the literature. 

In figure (^), equation ( 4.19| ) is plotted versus pq, for various positive values of bg. In 
all cases, di is an oscillatory function of po, and there are multiple regions of scalarization. 
We have verified numerically that scalarization does actually occur when di < 0. However, 
application of of the stability criterion of §3^ shows that these scalarized stars are unstable 
whenever 5^ > 0. 
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Figure 2: M vs R in Brans-Dickc tlicory for various equations of state and values of Oj. The 
central value, po = Pq/pq varies along each curve. The starting point of each curve {po 0) is at 
M = R = 0, and the endpoint of each curve corresponds to the ultra-relativistic limit where po — ^ p™'^^- 
The values of (j)o are chosen such that (jfoo is constant along each curve. The stellar configurations 
become unstable after the first turning point where dM/dR = 0. Top: Incompressible stars, for which 
pmax _ ^ Notice that the curves with non-zero scalar coupling have smaller maximum radii and 
masses. Bottom: Relativistic polytrope models of neutron stars, as defined in [Q, for which p™^'^ is 
finite (colour online). 
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Figure 3: Top: A comparison of ^ = a^/a^ = Q/{Mas) vs s = GM/R for incompressible stars in 
Brans-Dicke theory, for various values of flg. The curves start at s = 0, ^ = 1 where po ^ 0, and 
terminate where pq oo. Bottom: The same quantity comparing incompressible stars with two kinds 
of neutron-star equations of state (relativistic polytropes) given in ||l8|, using Brans-Dicke theory to 
leading order in ag. 
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Figure 4: The coefBcient di plotted versus po, for constant-density stars, for several choices of bs < 0. 
Scalarization becomes possible for bg < —4.329. 




Po 

Figure 5: The coefficient di plotted versus po, for constant-density stars, for several choices of bg > 0. 
There are multiple regions of scalarization. 
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4.3 Non-relativistic polytropes 



As an example for which the scalar-field dependence of the above constraints can be more 
explicitly explored, consider the case of Newtonian polytropes discussed §3^ with equation 
of state P = Kp^+^/^. In this case for weak central coupling the scalar field, and the 
dimensionless density, 9 = (p/po)^^^, can be expanded in a^, bg, as well as the parameter 
C = (1 + x)Kpl^'' = (1 + x)Po/po = (1 + x)Po: 



i i,j i,j,k 



(4.21) 
(4.22) 



We show below that C ^ 1 for the polytropes of practical interest, such as white dwarfs and 
main-sequence stars. 

In the limit = 0, eqs. ( 3.73 ) - ( |3.74| ) become 



26' 

w 
2ip' 

w 



0^ 



(4.23) 
(4.24) 



where for later notational convenience we use primes in this section to denote differentiation 
with respect to w. Equation (4.23) is called the Lane-Emden equation, and its solutions are 
well-studied because of its important role in the theory of stellar structure. It can be solved 
analytically when x = 0,1,5 [^]. The solutions of the Lane-Emden equation with initial 
conditions ( 3.68|) are called Lane-Emden functions, and are denoted by Q^{w). Thus, 

^(o)H = 0xH • (4-25) 
If we now take 6^ = then equation ( [4.24 ) can be solved for ip in terms of 9, giving 

^'^^^^^ = _i_e-C(-©x)' y'"(u,20'j'ef(-ex)'d^ (4.26) 

= -- + -^e-'^^^^-y r e^(-®x)'du; . (4.27) 
w J 

Expanding this in powers of C then yields 

mo)^C(l-e^) + 0(a. (4.28) 

and so 



¥^(0,0,0) = and 92(0,0,1) = 1 " ©x 



(4.29) 
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Expanding the matching equations, eqs. (|3.79| ) - ( 3.81 ) in powers of Oq, bs and C similarly 
yields 



GM 

Q 

Mao 
- 00 



ao 



-CWQ\{W) + 0{al), 

c 



1 - Q^iW) - WQ\{W)\ + 0{al) + 0{hs) + O(C') 
Notice in particular that eq. ( [4.31 ) implies that 

aA^ao\l + 0{al) + 0{bs) + 0{C)\ , 



(4.30) 
(4.31) 
(4.32) 

(4.33) 



which is consistent with the limit a a — ^ o(0oo) for weakly self-gravitating stars. Similarly, 
equation (4.30) can be re-written as 



GM 



-CtsW^Q'^ + 0{al) 



W'Q'JW) 



^(3-x)/2x 

Po 

^1 



4ttG 



+ 0{al) . 



(4.34) 



which is independent of po, as advertised, when x = 3. 

To get an idea about the validity of the C expansion, we close this section by estimating 
its size for the examples of white dwarfs and main-sequence polytropes. 

White dwarfs 

A white dwarf can be modeled by a degenerate fermion gas which satisfies in the ultra- 
relativistic limit the polytropic equation of state, eq. ( p. 62 ), with 



K 



31/3^2/3 



rribj 



and X = 3 . 



(4.35) 



Here Yg is the number of electrons per nucleon and mj, is the average nucleon mass |61]. 
Although the fermion gas is highly relativistic, the gravitational field generated by it is not 
strong, allowing use of the non-relativistic expressions developed above. 

The typical central densities of white dwarfs are pq ~ (lO'^ . . . 10^"^) kgm-3 |6]], corre- 
sponding to 10~^ ^ C ^ 10~^ for Ye = ^. This shows that a perturbative expansion in ( is 
likely a good approximation for white dwarfs. 

Since the ultra-relativistic limit of the fermion gas is used, the mass calculated below 
using the matching conditions — i.e. equation ( |4.34| ) — is greater than the actual mass of 
the white dwarf. It is instead to be regarded as the Chandrasekhar limit: an upper bound on 
the mass of white dwarfs. Because x = 3 its value turns out to be independent of po, and for 
Ye = ^, its value is approximately 1.4Mq [32|. 



- 28 - 



Main-sequence models 



In the Eddington stellar model, a star is regarded as an ideal gas whose energy is transported 
by radiation, and it is assumed that the gas makes up a fixed fraction of the total pressure, 
(3 := Pgas/P = const. Here Pgas is the pressure of the ideal gas, and the total pressure is 
P = Pgas + -Prad) where Prad is the radiation pressure. The Eddington model leads to a 
polytropic equation of state with 



K 



45 /P 



1-/3 

/34 



1/3 



and X = 3 , 



(4.36) 



where Rg is the universal gas constant, fi is the molar mass of the ideal gas and /c^ is Boltz- 
mann's constant [^]. Main-sequence stars can be approximately described by the Eddington 
standard model, even though convection also plays a role in heat transfer for more realistic 
models. 

For a more accurate single-polytrope model of the Sun, x = 3.35 and C ~ 10^^ |6^], so 
the perturbative expansion in remains a good approximation. 



5. Incompressible stars 

We now specialize the discussions of the previous sections to the special case of an incom- 
pressible star, for which the stellar density, p, is constant. Since constant density can only be 
consistent with the pressure gradients required for hydrostatic equilibrium if p 7^ Pip)y we no 
longer impose this kind of equation of state. (It is not necessary in any case, since the closure 
of the field equations is now accomplished by the incompressibility condition, p{r) = pQ.) 

The purpose of this exercise is to have a toy example for which all of the above ma- 
nipulations can be simply carried through explicitly in closed form. Performing the same 
exercise for GR provides an interesting example that displays the main features of relativistic 
structure, including the existence of a maximum compactness for a star, S(o) < |, that can 
be supported against gravitational collapse. The maximum that is found for incompressible 
stars turns out to provide an upper bound to the maximum compactness that can be achieved 
with other equations of state. 

Several earlier works have numerically investigated incompressible stars in scalar-tensor 



gravity 40, 45|, but our quasi-analytical treatment of these stars is new. 



5.1 Incompressible stars with quasi- Br ans/Dicke scalars 

We cut right to the chase and specialize directly to qBD scalars, for which a{cf)) = + 6s</>, 
since this case is broad enough to be of wide interest, but restricted enough to be explored 
in detail. 
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Field equations 

Taking the equation of state to he p = po, or q{p) = 1, equations ( p. 39 ) - ( 3.41| ) become 



At 
P 



.A + i ^iaMi+bs^m + al u{l - 2^)02 



2u{l - 2^) 4(1 - 2^) 



(3 - 4^)9i e^al^{l+bsV/2) 



+ 



2ii(l - 2;u) 8u(l - 2fi) 
Similarly, the function f{p) defined in eq. ( p. 15 ) becomes 

l+P 



(1 + bsip)il - 3p) + 2(1 - p)u^ 



(5.1) 
(5.2) 

(5.3) 



f{p) = In 



1+Po 



(5.4) 



and so the baryon number density, n(r), computed from eq. ( |3.33| ) is also constant, n = uq. 



The function Ai defined in eq. (3.43) similarly becomes 

fU 



M 



du. / — - — 
V 1 - 2m 



exp 



Q 2 I t , 

3ao¥' ( 1 + — 



(5.5) 



5.2 Perturbative solutions: leading order 

Because the GR problem can be explicitly integrated for incompressible stars, equations ( |5.1D 
|) can be solved analytically when Oq = 0. The zeroth-order profiles, ^(q) ^^'^ P{o)) 



given by |65, 36 1 



u 



(l + 3po)Vl- V3-(l+po) 
3(l+po)-(l + 3po)Vl-V3' 



(5.6) 
(5.7) 



as a function of the central pressure pq. We return to computing the profile, (p{u), below. 

At leading order the stellar radius is determined as the zero of p(o) > which vanishes at 
u = ?7(o) where 

_ 12po(l + 2po) 
(l + 3po)2 ' 

corresponding to the radius r = i?(o) with 



(5.8) 



(0) 



/3po(l + 2po) 



A^(l + 3po)V ^^'^/'o 



(5. 
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The leading components of the functions relevant to matching to the exterior solutions — i. e. 
T ^ s and M. — are 



-^(0) = ^.(0) - 12 (3^) In (i^) ^.(0) , (5.10) 
_ 2po(l + 2po) 

(l + 3po)2 ' ^^-^^^ 

A^(0) = 3V3 jarccos J - 2 1 , (5.13) 

where, as before, the subscript denotes evaluation at u = U, so 93^,(0) := V'(o)(f^(o))- 

Equation ( 5.12| ) implies that the compactness is an increasing function of po, which 

vanishes when po = and asymptotes to | as — ^ c«. Thus we reproduce the GR result 

< S(o) < I for constant-density stars. This prediction gets modified once ©(oq) corrections 

are included, however, as is discussed in detail in the next section. 
The mass-radius constraint, eq. ([4. 81), in this case becomes 



s^^^'-^AlR' + 0(al), (5.14) 

which states M oc R^, as might be expected for constant density. This behaviour is seen 
explicitly in the = curve in Fig. (^). 

The scalar-coupling constraint, eq. ( [4.11 ), similarly becomes 

«A = -, tPt^ ^ + 0{al) , (5.15) 

1 + bsif^o) - 66s (1 - 2s) ln(l - 2s)(^^(o) 



where 



Po= ' ';,_'^ + g(«o) (5.16) 



should be substituted into (/7*(o) and 9?*(o) on the right-hand side. Since further progress 
requires knowing the scalar profile, we next turn to solving its field equation. 

Scalar profile in Brans-Dicke theory (vi^hen bs = 0) 

If bs = 0, then the model reduces to Brans-Dicke theory, and a((/)) = oq = a^. Equation 
( |5.3D with Oq = becomes a first-order linear differential equation for ip, which can be solved 
analytically: 

27(1 +po) f '"";^ - a/T^Vs) - 4(1 + 3po)u 

ifm) (u) = , ^ ^ ^ == . (5.17) 

' 12n7r^V3(3(l+po)-(l + 3po)x/r^V3) 
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Integrating this expression once more with respect to u then gives 
1 



8(2 + 3po) 



(9po + 5)(3po - l)ln (3(1 +po) - (1 + 3po)V'^-u/3 

arcsin 



9(l+po) l + 3po + 3(l+po)\/l- V3 



9(l + 3po)2(l+Po) 



Li2(A_) -Li2(A+) + ihi 



1-A, 
1^ 



where 



and 



16(2 + 3po)y^ 

a + 3po)(Vl-n/3 + iyV3) 



\A73 

arcsin u/3 



A± : 



3(l+po)±2V2 + 3po 



Li2(z) = 



00 1. 



k=l 



, (5.18) 

(5.19) 
(5.20) 



is the dilogarithm function. 
By using the identity 



psmx 



k=l 



1 — 2p cos X + ' 



(5.21) 



it can be shown that the imaginary part of ( 5.18| ) is constant, and can thus be absorbed into 
an integration constant. The final normahzed and manifestly-real expression for (/?(o) is then 



1 



8(2 + 3po) 



(9po + 5)(3po - 1) log Q(l +Po) - ^(1 + 3po)Vl- 

arcsin y^u/3 



9(l+po) l + 3po + 3(l+po)\/l- V3 



9(l + 3po)2(l+Po) 
16(2 + 3po)3/2 



arctan 



2^(2 + 3po)u/3 



l + 3po-3{l+po)^l-u/3 
+ 5}?[Li2(A_) -Li2(|A_|) -Li2(A+) +Li2(|A 



u/3 

arcsin a/ u/3 



+ 4(1+^0) 



(5.22) 



value 



If »o > ^) then the inverse tangent in the above expression changes branch at the critical 

v3 



4(2 + 3po) 

^crit = TTTT- T2 < '^(O) ; 



3(1 +Po) 

SO that (^(0) is a continuous function of u. 

Substituting equation ( |5.17] ) into ( [5.11 ) yields 

9(l+po)(l + 3po)' 

^(0) - 



(5.23) 



16[po{l + 2poW/^ 



■ arccos 



l+Po\ Alpl + 34po + 9 



l + 3poy 8po(l + 2po) 



(5.24) 
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so the scalar-coupling constraint, eq. ( 5.15 ), can be explicitly evaluated, 
'5 9 9\/r^^[27s-14-9^/r^^(2-3s)] 



ao 



2 4s 
6 



8[s(5 - 9s + 3V1 - 2s)]3/2 



arccos 



VI - 2s 



ao[l-^s + Ois^)]+0{al) 



(5.25) 
(5.26) 



This confirms that although — >• oq in the non-relativistic limit s — t- — consistent with 
eq. (|3.37|) — it is in general depressed relative to ao for relativistic systems, even when 
additional powers of are neglected. 

Notice that in the opposite limit we have 



lim a A 



^^/2arccos(l/3)-^^ +0(c 
~ -0.35900 + 0(0;]). 



(5.27) 



Equation ( |5.25| ) is plotted in Figure ^, showing that for incompressible stars passes through 
zero, changing sign at s ~ 0.398. This is compared in the same figure to the corresponding 
curves for neutron stars modeled by relativistic polytropes, using the equations of state EOS 



A and EOS II defined in |18]. For small s all three curves agree reasonably well, but differ 
for large s. For large s the neutron star curves significantly deviate from the constant- 
density curve, as might be expected given that relativistic polytropes have a maximum value 
^max _ ^ _ 2^ where 7 is the polytropic index, while constant-density stars have no such 
maximum value for pQ. 

Substituting equations ( ^.22 ) and ( 5.17 ) into (|5.10|) yields 



(0) 



1(1 +P0 



41pg + 34po + 9 , 

- 8..n^w + 3^°) 



+ 



8po(l + 2po) 

9(l+po)(l + 3po + 3pg) f 1+po ^ ( l+Po 

= — = log(l + po) — arccos 

4(2 + 3po)VPo(l + 2po) ^^^50(1 + 2^0) Vl + 3po, 



9(l + 3po)'(l+Po) [ fl+Po" 

: — s arccos 

16(2 + 3po)=^/2 \ \l + 3po, 



' 2 + 3po ^ 
^Po(l + 2po). 



3/2 



log 



1+Po 
.1 + 3po. 



/ 2^Po(l + 2po)(2 + 3po) 

+ arctan 7^ 

' 1 - 3p2 



where 



- K[Li2(A_) - Li2(|A_|) - Li2(A+) + LiadAn 



^ ^ l+Po + 2Vpo(l + 2po) 
^ 3(l+po)±2V2 + 3po 



(5.^ 



(5.29) 
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is the value of X± when u = U(q), and 



arctan X 



(5.30) 



(0)' 

{arctan X — ir ifO<po<-^i 
arctan X if po > • 

Equation ( ^.28 ) is plotted in Figure and is compared to the corresponding curves for 
neutron stars. Again, the curves are close for small s, and diverge for large s. J^(o) (s) is positive 
and increasing on the interval < s < 0.25; positive and decreasing when 0.25 < s < 0.36; 
and negative and decreasing for 0.36 < s < |. The maximum value, attained at s ~ 0.25, is 



-^(0),max ~ 0.27. As S 



|, J^(o) tends to — oo. 



Solution for 6^ / 

Next suppose bs / 0. To calculate i^(o), change variables from ((^(0),^) to {ip,z) with 



= 1 + 6s95(o) > ^ 
Initial conditions for ip{z) are then 



1-^ 
3 



6^(1 - 3po) 



(5.31) 



(5.32) 



2 = 
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Equation (5^) with Oq = becomes 

d^ip / 7 5 e \ dip aj3 z — q 



where 



dz^ ~^\z~^z— l^z — aj dz ~^ z{z — l){z — a)'^ ' ^ ^ 

7 = ^ = 1, e = 1 , (5.34) 



l + 3po' 2 V3Po + iy ' ^ ' 2' 



« = l > /3 = l l + , (5.35) 



and 



7 + ^ + e = a + /3 + l. (5.36) 

Equation ( |5.33| ) is called Heun's equation [^], and is a linear second-order differential 
equation with singularities at z = 0, l,a,oo. It is a natural generalization of the hyperge- 
ometric equation to the situation having four regular singular points. The solution which 
satisfies initial conditions ( 5.32| ) is the local Frobenius solution about z = with exponent 0, 



oo 



and is given by the power series 

V'(z) = HeunG(a, q\ ot, (3, 7, 5\ z) = c^z^" , (5.37) 

where the first two coefficients are given by 

Co = 1 , ci = 4:: , (5.38) 

07 

and the higher coefficients are found by solving the three-term recurrence relation 

(r-l + a)(r-l + /3) Cr-i 

-[r(r - 1 + 7)(1 + a) + r{a5 + e) + ^] Cr 

+a(r + l)(r + 7) c,+i = 0. (5.39) 

In terms of hg and po) the recursion relation for the coefficients Cr become 

co = l, ci = 6,(1 - 3po) , (5.40) 

2(1 + 3po)(?'^ + r - 2 + 66,) c^_i 

- [r(6po?" + 9po - 1) + 36s(3po - 1)] 

-(r + l)(2r + 3) c^+i = 0. (5.41) 

This implies the coefficients Cr of the power series ( ^.37]) can be written explicitly as a poly- 
nomial of degree r in 6,, 

r 

Cr=Y.<^^s, (5.42) 

i=0 
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(r) 

where a- is itself a polynomial in po of degree r. 

The solutions for a[ and Cr are found explicitly in the Appendix, where it is also shown 
that the coefficient of the largest power of bg has a particularly simple form: 



[6(1 - SpoW 
(2r + 1)! 



(5.43) 



Because 6s is relatively poorly constrained, it can be larger than unity so far as phenomenology 
is concerned. In this case eq. ( ^.431) can be used to obtain an approximation for ip{r) for large 
bs- This gives (see Appendix for details) 



A 



(0) 



1+Po 



Po 



1 

261 



^cosh vt - '-^^] + y y p^^}p^)fw{T) 



fc=i i=o 



(0) 



2po 



k=l j=0 



sinh VT l+po 
2^ 



/T 
Po 



L cosh Vt — 1 



Tfk,j+i{T) 



Pk,j{po) 



Lk+l 



[6(1 - 3po)]2'= ' 



(5.44) 



(5.45) 



where T = 66sPo(l — 3po)/(l + 3po) and L = log(l + po) — log(l + 3po)- If bg > 0, then 

T < (6 - 4^/2)6s ~ 0.346s. If 6^ < 0, then T > -(6 - 4^2)16^1 0.34|6s|. The above 

expressions (|1D and (|5l5|) can be used to calculate the second constraint ( [4.11|) . 



Compactness vs central density 



Equation ( ^.12 ) describes how the compactness depends on po in GR. In order to find how 
scalar-matter couplings modify this behaviour, it is necessary to solve the equations of stellar 
structure to first order in Cq, and calculate S(i). 

Figure ^ plots the compactness vs po in Brans-Dicke theory, for various values of the 
Brans-Dicke coupling gq = a^. Notice that the compactness eventually stops growing with 
Po, approaching instead an asymptotic value as po — ^ oo. In GR, this asymptotic value is 
GM/R = 4/9, which is the maximum allowed by Buchdahl's theorem. As oq increases, this 



asymptotic value decreases. This is consistent with the results of [l40|. 



5.3 Perturbative solutions: next-to-leading corrections 

In this section, O(ao) corrections are calculated. The defining equation for U is p{U) = 0. 
Expanding it in powers of Oq yields 



(1) 



12(1 +Po) 
(1 + 3po)^ 



P{i)i.U(^o)) ■ 



(5.46) 
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Figure 7: The compactness s = GM/R plotted versus po = Po/po for constant-density stars in 
Brans-Dicke theory, for various values of Oq. 



Solving equations (5J) - ( |5.3| ) perturbatively in Og yields 



pi(u) 



1 

/ dnV^ 
/n Jo 



f(0) 



n(l - u/3)(vi(o))' 



m + 2 

2(1 +po)v9(o) (1 + bsV?(o)/2) 

(l + 3po)Vl- V3-3(l+po) 
2(l+po)/i(i)(n) _ (1 + po)(3 - 2n) yi - ^3 - (1 + 3po) 

- ^x/3 " ((l + 3po)Vl- V3-3(l+po))2 
2(l+Po)Vl- V3 r du 

6((1 + 3po)x/r^W3 - 3(1 +po))2 70 (1 - n/3)3/2 



(5.47) 



+ 



J(n) 



(5.48) 



where the function J{u) appearing in equation (|5.48 ) is given by 

J{u) = 2 [6(1 + 3po)\/l- V3- (1 + Po)(2u^ -9u + 18)] u(l - u/3)(</i(o))^ 

^s93(o)" 



-(l+Po)(4u"-18tx + 9)(^(o) (^1 + 
The perturbation to the scalar profile is similarly 



(5.49) 



(5.50) 
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where the functions and ^|) are local Frobenius solutions of equation ( ^.33 ) (with parameters 



( |5.34| ) - ( 5.35D ) about z = with exponents and — ^, respectively. They are given by 



ih = HeunGi , -bg ■ : 



3(1 _ V1-8V3), ^(1 + Vl-86s/3), y ; , (5.51) 



^ = ' Heuncf 6U^p, - 1) + 1- 6p, 

Vl + 3po 4(l + 3po) 



1 + ^ Vl - 86s/3, 1 - ^ Vl - 86./3, ^, ^; , (5.52) 



where z = (1 — -y/l — u/?>)/2. The coefficients <I>i and <I>2 are given by 



$1 = 288 / + +^)(1 -2z)2[z(l -z)]3/2[i^ (1^3^^)^]^^^ (5^53) 

Jo 

$2 = -288 / il'iFp^i) + + H){1 - 2zf[z{l - z)fl^[l + (1 + 3po)z] , (5.54) 



where 



G 



H 



1 + 6^^(0) + 8z(l - z)93(o) 
32z(l - 2)(1 - 2z)2 
(4(1 + 3po)z + 1 - 3po)(l + fe,y9(o)) 
48z(l - z){l - 2zf{l + (1 + 3po)^) 

[12(1 + 3po)z3 - 6(1 + 7po)z^ + (9po - 5)z + l]^(o) 
12z(l - z){l - 2z)4(l + (1 + 3po)^) 
(4(1 + 3po)z + 1 - 3po)v9(o)(l + bsV?(o))(l + bs>^(o)/2) 
2iz{l -z)il- 2z)2(l + (1 + 3po)z) 



(5.55) 



(5.56) 



(2(1 + 3po)z + 1 - po)v?(o) (1 + 6^9^(0) /2)93(o) 

(l-2z)2(l + (l + 3po)^) • ^ ' ' 
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Figure 8: Comparison of /i vs u calculated perturbatively and numerically for an incompressible star 
in quasi-Brans/Dicke theory with Cq = 0.1 and 6s = 4, po = 0.1 (top left); 6^ 4, po = 1 (top right); 
hs = —4, po = 0.1 (bottom left); and hs — —4, po = 1 (bottom right). All curves terminate at the 
stellar exterior, u ~IJ . 



The first-order corrections to the external parameters are given by 

3(1 +Po)' 



M 



^{1) = '/'(I) - 12e2^L0(i) 



+12 1 + 



2po(l + 2po) 
(l + 3po)' 



2,0(1 + 2,0)^^' ^'^'^^^^^^W 



2lPo( 1 + 2po) 
+ 3po)2 



2po(l + 2po) 



L 



A io2L- , 3(1 +po)^ , ^ .(l + 3po)' 



2po(l + 2po) 
-1728e^^^gii±^(vi(o))^ 



po(l + 2po) 



95(0)^(1) 



(1 + 3po)2 

5(1) = ^(1) + 2e^^p(i) - 144L£^ ^ ^^^^^ 



IS 2L + 4po + 1) . 
18e ^^V5(o)P{i) 



Po(l + 2po) 



2L„ ....Po(l + 2po)^4L(^^^^)2^ 



24(l+po)v/3po(l + 2po) , 12V3po(l + 2po) 

« = OT^^^? ''^^^ ^ ^0 ^^'^ 

duj- + hs^(o)/2) + 2n(l - n/3)(9i(o); 



(5.58) 

(5.59) 
(5.60) 



10 V 1 - V3 

where L = log(l +po) — log(l + 3^0), and the profiles are all to be evaluated at t/{o)- 



(5.61) 
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Figure 9: Comparison of p vs m calculated perturbatively and numerically for an incompressible star 
in quasi-Brans/Dicke theory with Oq = 0.1 and bs — 4, po = 0.1 (top left); 6^ = 4, po = 1 (top right); 
bs = —4, po — 0.1 (bottom left); and bs — —4, po = 1 (bottom right). All curves terminate at the 
stellar exterior, u — U. 



If bs 7^ 0, then the relation a a = dlogM/d(f)oo can be used to simplify Ai 



(!)• 



, . 9^3 f l+Po\ V3po(l + 2po)(41pg + 34po + 

= -r^ — arccos ' > ^ 



2bs \l + 3poJ 6,(l+po)(l + 3po)2 



96(l+po)V3j^g(l + 2po)3 . , , , , 8V3po(l + 2po) 

6.(1 + 3po)^ ^(°) + + TT^o '''' 

, 24(l+po)V3po(l + 2po) 

+ (TT^^ ^'-''^ 



5.4 Comparing perturbative solutions with numerics 

Part of the utility of analyzing the incompressible star in such detail is that such explicit 
expressions for the perturbative solutions allow a detailed comparison with direct numerical 
integrations. This helps indicate the domain of validity of the perturbative results. 

First, we look at the profiles for the physical variables n{u), p{u), and ^p{u) across the 
interior of the star. The quantities fi, p and ip are respectively plotted versus u in Figures 
, for Oq = 0.1 and several choices for bg, and pQ. The line labelled "Zeroth 
Order" plots the zeroeth-order result, e.g. /i(o), while the line labelled "First Order" includes 
also the first correction, e.g. //(g) +'^oA^(i)- Notice that the curves all lie close to one other for 
small u, but begin to separate at the stellar exterior, u ^ U. Furthermore, the separation is 



(|), (|) and (10 
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Figure 10: Comparison of tp ys u calculated perturbatively and numerically for an incompressible 
star in quasi-Brans/Dicke theory with Oq = 0.1 and bg = 4, po = 0.1 (top left); bg = A, po = 1 (top 
right); 6s = —4, po = 0.1 (bottom left); and bs — —4, po = 1 (bottom right). All curves terminate at 
the stellar exterior, u = U. 



largest for the more relativistic stars, for which po is larger. However in all cases displayed 
the perturbative results capture the full numerics quite well throughout the entire star, with 
the strongest deviations happening for ip{u) when 6^ < 0. 

Of more practical interest is a similar comparison of the accuracy of the perturbative 
expressions for plots that directly relate observable quantities to one another, such as plots 
of a A vs s. Examples of these are given in figures (11), which give s, A = 0^/00, T = 
{(poo — 4'o)/cLo and M. as functions of the central pressure, po = Pq/po for the special case 
of Brans-Dicke theory {bg = 0) with Oq = = 0.1. These again show good agreement 
between perturbative and numerical calculations, with the biggest deviations arising in the 
most relativistic settings (largest po). 



6. Conclusions 



In this paper we set up the equations of stellar structure, with the stellar interior modeled 
as a spherically symmetric, static fluid, and with gravity described by a scalar-tensor theory 
with a single light scalar coupling to matter only through its coupling to a Jordan frame 
metric. For practical reasons, and for the purposes of making contact with earlier workers, 
we focus on the special case where the scalar-matter coupling function does not vary strongly 
with the field, a{<j)) ~ + bs4>- 
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Figure 11: Comparison of pliysical quantities as functions of central pressure, po = Po/po, calculated 
perturbatively and numerically for an incompressible star in Brans-Dicke theory (bg = 0) with = 0.1. 
The plots show compactness, s (top left); external coupling, A = aAja-s = Q/^Ia-si (top right); 
^ = (0OC — 4'o)/o,s (bottom left); and M (bottom right). 

We seek solutions to these equations, for a variety of equations of state, in the special 
case where the scalar-matter coupling at the stellar center is small, Oq = a^((/>o) <C 1. We 
obtain solutions as perturbations to those of General Relativity. By comparing these solutions 
with explicit numerical integration we verify that the perturbative approximation works well 
throughout most of the star. 

These perturbative solutions have the merit of being very simple to integrate numerically, 
and of allowing analytic solutions for some choices of equation of state. This is very convenient 
for efficiently exploring different choices for the scalar properties, and scalar-matter couplings. 

We use these solutions to compute the form of the observable relations that are imposed 
among the external properties of the stars by the condition that they match continuously to 
the stellar interior. There are two such relations among the four external variables, M, R, 
Q and and our semi-analytic approach allows a simple exposition of how these relations 
depend on scalar properties. These properties ultimately underly any tests of scalar-tensor 
theories using astrophysical systems, such as binary pulsars. 

Finally, these methods are applied to the illustrative case of an incompressible star, for 
which the density is constant. In this case the solutions generated by the oq expansion may be 
found analytically, making the comparisons with numerical results particularly simple. Again 
we find that the perturbative expressions agree well with the solutions obtained by numerical 
integration. 



- 42 - 



Acknowledgements 



We thank Nemanja Kaloper and Maxim Pospelov for useful discussions. This research was 
supported in part by funds from the Natural Sciences and Engineering Research Council 
(NSERC) of Canada. Research at the Perimeter Institute is supported in part by the Gov- 
ernment of Canada through Industry Canada, and by the Province of Ontario through the 
Ministry of Research and Information (MRI). 



A. Expansions for large bs 

This appendix evaluates the large-6s expansion for the properties of incompressible stars. 



The coefficients Cr of the power series (5.37) can be written explicitly as 



Cr = ^a^^K, (A.l) 



i=0 

where af"^ is a polynomial in po of degree r. It follows from equation ( |5.41| ) that a^^ satisfy 
a two-term recurrence relation, which can be solved explicitly, giving 

"'^ " (2r + l)! • ^ > 

(r) 

The coefficients Oq vanish for r > 1, and the remaining coefficients are given by 



where Pkj is a polynomial of degree 2k with rational coefficients. The first of these polyno- 
mials are given by 

^^1,0(^0) = 2(63pg - 18po + 7) , (A.4) 

^'i,i(Po) = 8(9pg + 9po + 4), (A.5) 

-P2,o(Po) = |(-3321p^ + 2916pi^ + 6426p^ + 4644po + 679) , (A.6) 

P2AP0) = |(5913p^ - 21708p{] + 29502p2 + 17028po + 2593) , (A.7) 

^^2,2(^0) = f (729p^ - 729pl - 1989p2 - l251po - 176) , (A.8) 

P2,3{Po) = 32(9p^ + 9po + 4)2 , (A.9) 
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and the higher ones can be computed from the relation 

, 21-1 



i=0 

^ g (2, + l)0 + 2)(iy(l-3poft-| g _ 

j=l 1=0 

^ (1 - 9po - 6jpo)(l - 3po)jj! D . . 4 
+ 2^ y—jT, 2^ ^'«-i,.(po)j 

+4 E (^^+;Hl + 3po).! gp,-.,fa)(i-l)-, (A.10) 



which holds for Z > 3 and r >l + 1. 
From the identities [|67|1 



^(-l)^' {k + af = , (n > 0, a G M) (A.ll) 

N 

^)(A: + a)"-^ = 0, (iV>n>l, aGM) (A.12) 



k=0 



it follows that 



Write 



[8(9pg + 9po + 4)]^ 

-Pa:,2A:-i(Po) = • [k > I) (A. 13) 

2fc 

PkAPo) = T.^klPo- (A.14) 

By means of telescopic summation, the identity ( A.10| ) can be used to derive recurrence 
;ions for the coefficients q^^j. In general, 
but for I = 2k, they are particularly simple: 



relations for the coefficients q^^j. In general, these recurrence relations are very complicated 



2k-l 

+ 2),i-) + ^ (z + l),if = 324,f , (A.16) 

1=2 

i=i + l V / 



where k > 2 and 2 < j < 2k - 2. 
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Now, (/5(Q) can be expanded in the form 



'/'(o)(^) = TT 



1 



6^ ^ 6? +••• • 



It follows from equations ( A.2 ) and ( [A. 3D that 



(A.18) 



where 



where 



sinh \/i 



1 



2fc-l 



[6(1 - 3po)] 



r=A:+l 

d 
It 



2k 



j=0 



r! • fr-' 



(2r + l)!(r-A;-l)! 

jk+iiiVi) , 



2 y 




2j. "+ 2 ^ ' \ ' \x dx I X 



1 d\^ sinx 



(A.19) 
(A.20) 



(A.21) 



(A.22) 



are the spherical Bessel functions, and t = 66^(1 — 'ipo)z = 36<j(l — 3po)(l — y^l — ti/3). 
The spherical Bessel functions satisfy the recursion relations 



xji-i{x) + xji+i{x) = {21 + IjjKa;) , 
lji~i{x) -{1 + l)ji+i{x) = {21 + l)j[{x) . 

Therefore, the fkj satisfy the recursion relations 

fk,i = fk+ifl + {k + l)/fc,o , 
tfk-1,0 - 4/fc+i,o = 2(2/;; + 3)fkfl ■ 



(A.23) 



(A.24) 



The first of these relations can be used to write fkj explicitly in terms of fi^o: 
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m=l 



j—m+1 



fk,j = fk+j,Q+Y] V, , {k + li)...{k + lm) 



fk+j—m,0 ' 



(A.25) 



for j > 1. The tilde on the sum means that ah numerical factors which come from overcounting 
should be deleted. For example, 



V, , ik + h){k + h) = {k + lf + {k + l){k + 2) + {k + 2) 
^ — '11,12=1 

Equations ( ^.lOD - dsTTll ) then become 



A 



(0) 



(0) 



1±^^ cosh^/T- 



2po 



00 2k-l 
k=l jr=0 



sinh vT 



00 2fc-l 



Pk,j{Po)fk,j+i{T) 




(A.26) 



(A.27) 



6^^6(1 - 3po)]2'= 



(A.28) 



where T = 66sPo(l ~ 3po)/(l + 3po) is the value of the variable t corresponding to u = U, 
and L = log(l +po) - log(l + 3po). If 6^ > 0, then T < (6 - 4\/2)6s ~ 0.346^. If bs < 0, then 
T > — (6 — 4\/2)|6s| ~ — 0.34|6s|. The above expressions (^^ and (lA^ ) can be used to 
calculate the second constraint ( [4.11|) . 

These perturbative solutions have a limited regime of applicability. They break down 
as bs — ?■ 0, because the l/bg expansion fails. But they also break down at the onset of 
spontaneous scalarization. 
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Table 1: Table of Notation 
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